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These examples show that the method is applicable to higher 
frequency ranges. 


NOMENCLATURE 


A 

A («) 


c" 

e r 

F 

f 

4 2) 


total dimensionless area of finite element domain 
area of element e 

modal amplitude of wave moving in +r direction away 
from origin 

modal amplitude of wave moving in — r direction 
towards origin 

dimensionless circular cylinder radius a^/L^ 
dimensional radius of circular cylinder (see superscripts) 
adiabatic speed of sound 
error, Eq. (26) 

column vector, right hand side of global matrix equation 
dimensionless frequency, Eq. (2) 

Hankel function of the first kind 
Hankel function of the second kind 

/=■ r 


J Bessel function of the first kind 

m 

K global matrix, Eq. (23) 

k wave number, Eq. (4) 

characteristic distance 

M number of elements in finite element domain 

M c number of modes needed for convergence to error e r 

M coe f number of modal coefficients used in eigenfunction 
expansion 

Mp tH number of grid points on interface S used in integration 
m mode number 

m mode number, Eq. (10) 

N number of nodes in finite element domain 

N^ local linear triangular interpolation function, N^(x,y) 

NS e >(x j)yj ) = (i = 1,2,3; j = 1,2,3) 

N S number of grid points on interface S 

n normal vector 

n outward unit normal vector 


t dimensionless time, 

W| global weight function associated with node i; 

Wi(x jiyj ) = S V] (i = 1 . . . N; j = 1 . . . N) 

W[ e ) local weight function associated with node i 
W m * interface weight function, Eq. (10) 
x dimensionless axial distance, x^/L^ 

y dimensionless transverse distance, y^/L^ 

ft angle between element outward normal and x axis 

6jj Kronecker delta (<5y = 1 for i — j; = 0 for i # j) 

e dimensionless complex acoustic permittivity 

$ angle between radius vector and x axis 

p dimensionless complex acoustic permeability 

p dimensionless perturbation (acoustic) density 

/»(*. y) = P*/Po,n* 

p 0 dimensionless ambient density, ^ 0 ^(x,y)/p or( , f # 

(p 0 (x,y) — 1 for uniform properties) 

ij> column solution vector, Eqs. (23) and (24) 

uft angular velocity, 2irf^ 

<jj dimensionless angular velocity, Eq. (2) 

Subscripts: 

a analytical solution 

I node location on interface S 

i node location in finite element domain 

0 mean fluid condition 

ref arbitrary reference location for normalizing material 

properties 

x x component 

y y component 

Superscripts: 

~ approximate solution 

# dimensional quantity 

' derivative with respect to argument 

J average value 

(e) element value 

1 incident component 


P dimensionless perturbation acoustic pressure, 

P # (x,y,t)/p 0 ref # c ref #2 

p dimensionless perturbation acoustic pressure 

p(x,y) = PCx.y.tJ/e -1 " 4 

Rj global residual error at node i 

r dimensionless radial coordinate 

S line interface between finite element and homogeneous 

domains 

S n partial sum of series with n terms 

S + region exterior to S 

S~ region interior to S 

s arc length parameter on S 


S scattered component 

METHOD OF ANALYSIS 

The present study is concerned with computing the acoustic 
scattering by a symmetric two-dimensional body of an impinging 
plane wave traveling in the -j-x direction. The spatial domain is 
divided into two subdomains, the homogeneous domain and the finite 
element domain, as shown in Fig. 1. The finite element domain 
contains a nodal grid system that covers the region of complicated 
geometry and material property variations. Linear triangular ele- 
ments are used and the subdomain interface is approximated by 
piecewise linear segments. In the finite element domain, an approx- 
imate solution for the total (incident + scattered) acoustic pressure 
p is calculated by the Galerkin method. In the homogeneous domain, 
which extends to infinity, an analytic solution (an eigenfunction 
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Figure 1 . — Finite element grid system. 


expansion) for the total acoustic pressure p a is derived by separa- 
tion of variables. 

The modal element method couples the two solutions by 
imposing continuity on the pressure and velocity at the interface 
between the two subdomains. This coupling results in a single 
matrix equation in which the eigenfunction coefficients and the 
pressure at the finite element nodes are calculated simultaneously, 
yielding a global representation of the acoustic Field. Generally, 
however, the desired quantity is the pressure field in the homoge- 
neous domain, so the main goal is to find the eigenfunction coeffi- 
cients; the finite element solution acts simply as a tool to obtain 
those coefficients. 

GEOMETRICAL MODEL 

A general finite element grid (Fig. 1) has been developed to 
study acoustic scattering from symmetric two-dimensional bodies 
with either hard (reflecting, rigid) or soft (absorbing) internal 
structures. The grid divides the finite element domain into rec- 
tangular, transition, and polar regions. The interior rectangular grid 
is used to model complicated geometries by a set of stepwise curve- 
fitting parameters, exemplified by Insert A (an aircraft inlet nacelle) 
in Fig. 1. The outer polar region is used to model more circular 
structures, such as an absorber-coated hull, as shown in Insert B in 
Fig. 1. The transition region links the polar and rectangular grid 
regions. The finite element domain can be restricted to the scatter- 
ing body, or can extend as far beyond the body as desired. 

For absorbing bodies, the entire body is covered with finite ele- 
ments to track the propagation of an incident wave penetrating it. 

It is generally accepted that 12 nodes per wavelength are required 
when using linear elements, so the higher the frequency of the inci- 
dent wave, the finer the grid must be. 

For rigid bodies with simple shapes, it is possible to reduce the 
finite element structure inside the body to a single ring of elements 
following the body’s contour, as shown in Fig. 2. This variation, the 


modal ring method, significantly reduces computation time and array 
sizes by eliminating nodes inside the body where the field does not 
penetrate. Coated rigid bodies could employ a small number of rings 
to calculate the pressure field in the penetrable coating. 

GOVERNING EQUATIONS 

Acoustic propagation in two-dimensional space can be modeled 
by the continuity, momentum, and state linearized gas dynamic 
equations in the absence of flow. For harmonic pressure propagation 
in a heterogeneous bulk material, Baumeisfcer and Dahl (1987) used 
these equations to derive the following dimensionless wave equation: 


_a 

dx 


1 dp 

d 

+ 

f 

1 dp 


dy 

K e dy J 


+ u 2 fip = 0 , 


( 1 ) 


where the dimensionless frequency associated with this Helmholtz- 
like equation is defined as: 


f = 


L*f # 

c ref# 


» 2*-f . 


(2) 


The notation is explicitly defined in the nomenclature. The refer- 
ence location is arbitrarily chosen; all nondimensional properties are 
normalized to their reference values. 

The relationship among acoustic “permittivity 11 £, acoustic 
“permeability” ft, and the physical properties of the medium is 
complicated in the general case. However, for propagation in 
homogeneous media without dissipation, such as air, e is the non- 
dimensional ambient bulk fluid density and ft is its inverse: 


* = /> 0 ( x >y) i 




1 

^ 0 ( x >y) 


( 3 ) 


For bulk absorbers, Baumeister and Dahl (1987, Eqs. (25) to (27)) 
employed Hersh’s model (1980) in explicitly relating e and ft to the 
porosity, a viscous loss coefficient, a heat transfer parameter, and an 
effective speed of sound of the medium. Briefly, for bulk absorbers, 
the real part of e is related to ambient bulk fluid density while its 
imaginary part is proportional to a viscous loss coefficient. The ft 


Homogeneous 



Figure 2. — Finite element ring grid system for rigid 
bodies. 
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term is inversely proportional to ambient bulk fluid density and the 
effective speed of sound, which depends on the heat transfer charac- 
teristics of the medium. The formulas for e and pi are given by 
Baumeister and Dahl (1987, Eqs. (5) and (7)). Morse and Ingard 
(1986, p. 253) also developed more general parameters for describing 
propagation in porous media for which e and pi can be related. 
For elastomeric materials, such as rubber, e and pi are generally 
complex because of the complex nature of the propagation speed in 
the medium (Capps, 1989, p. 10). In the present paper, the param- 
eters e and pi are treated as mathematical quantities independent 
of property correlations. Both quantities are complex with the 
imaginary parts associated with dissipative losses. 

In the homogeneous domain, the properties e and pi are 
assumed constant. Thus, by substituting Eq. (3) into Eq. (2), the 
governing equation simplifies to 


the +r direction and terms correspond to inward traveling 

waves. The scattered wave is then written as 

M coef -l 

p’ = £ A*H^ ) (kr)co 8 (m«) , ( 8 ) 

m~0 

where M coe f, the number of coefficients used in the expansion, must be 
set a priori. Hereafter, the eigenfunction terms are called modes, as 
commonly used in acoustic theory. Equation (8) satisfies the far field 
boundary condition — scattered waves move outwardly. The coeffi- 
cients are unknown, to be determined by the imposition of 

interface conditions. 

INTERFACE CONDITIONS 


^P a 



= o, 


(4) At the interface S between the finite element domain and the 

homogeneous domain, both pressure and velocity are continuous 
(Temkin, 1981, p. 80). The continuity of pressure 


where k 2 = w*fie. The subscript a is introduced to indicate that 
the acoustic pressure solution designated by p a is the analytical 
solution that applies in the homogeneous domain, in contrast to the 
pressure p designating the pressure in the finite element domain. 

ANALYTIC SOLUTION 

In the homogeneous domain, an exact eigenfunction expansion 
can be derived from Eq. (4) by separation of variables. This expan- 
sion, valid only in the homogeneous domain, is linked to the finite 
element solution, valid only in the finite element domain, by the 
modal element method. 

Consider an acoustic pressure plane wave traveling in the +x 
direction, striking a symmetric two-dimensional scatterer. In the 
homogeneous domain, the expression for the total pressure p a can 
be broken into the incident plane wave and the scattered wave: 


Pis’ - Pis- “ 0 


(9) 


can be expressed numerically by a collocation procedure (Lee and 
Cendes, 1987) or an integral weighting procedure (Baumeister, 1986b, 
Eq. (45)). The latter is used here with weight functions cos(m 0), 
so the continuity of pressure at the interface is expressed by 


9 - 2 jt 

fs [ p * - P ¥ S ■ 0 • 

9=0 

( 10 ) 

W . = cos(m’tf) 
ra ' ' 

(m* = 0 > l,2,...,M coef -l equations) 


where p a is the incident wave and p* is the scattered wave. A 
plane wave moving in the -{-x direction has the form 

= e ikx = e^co-W = j o(kr) + 2 £ i m cos(m0)J m (kr) , 

m=l 

when the time dependence is assumed to be 

The eigenfunction expansion in two dimensions is written in 
the following form (Morse and Ingard, 1968, p. 401): 

M co«f“ ^ . 

Pa “ E A rn H m ( kr )cos(m0) + A m H m^ ’(kr)cos(m0) . ( 7 ) 

m~0 

Equation (7) is valid only for problems that are symmetric with 
respect to the x-axis because of the assumed mathematical form of 
the incident wave in Eq. (6) (Eq. 1.2.9 Morse and Ingard) which 
leads to the cosine term in Eq. (7). Due to the form of the time 
dependence, terms correspond to waves traveling outward in 


In general, S can be parameterized by a function of r and 0, 
introducing a 0 dependence into ds. Unless S is a circle, then, 
the orthogonality of the cosine functions cannot be used, and a com- 
pletely numerical approach must be adopted to determine equations 
for the A^. 

For the problems considered so far, it suffices to apply a simple 
(6) quadrature to obtain acceptable results when approximating Eq. (10). 
S is divided into subintervals centered at points (rj,0j), which cor- 
respond to finite element nodes introduced later in the paper. Once 
M coef has been determined (based on the incident wave frequency and 
data in Figs. 3 and 4), the grid is set up so that the number of nodes 
on S is N g > 12M coef to adequately resolve the pressure field. 
Experimentation has determined that, at least for simple geometries, 
not all of these points need to be used in the quadrature to approx- 
imate Eq. ( 10 ) accurately. In the examples presented later, the 
number of points used is Mp tg = 2M coef . These nodes are equally 
spaced around a circular interface S, so every sixth point is used in 
the quadrature below. 

It is highly desirable to reduce the number of nodes used; doing 
so reduces the bandwidth of the solution matrix and yields a large 
savings in core storage and computation time. For more compli- 
cated problems, however, it will probably be necessary to use a 
better quadrature to approximate Eq. (10), which may increase com- 
puting requirements. 
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Figure 3. — Number of modes required to resolve the 
scattered acoustic pressure field on a cylinder of 
radius a subjected to a plane wave impingement for 
a relative error e r . 


Applying the midpoint rule to the chosen nodes gives 


Mpt. 


Es |P*( r I’ fl l) “ Pi] c os(m*0i) ri Atfj = 0 . 


I=l 


(m* = 0,1,2,...,^^— 1 equations) 

By expressing p a in terms of the modal coefficients in Eq. (8), 
Eq. (11) can be written explicitly as 


Mpu 

Es r l A0 l c °s(m*9j) 
1=1 


ikij cos 0j ^ * + TT (lJ/. x x 

e + Z A m H m ( kr l) COS Ml) “ Pi 




(i), 

m m 


( 12 ) 


- 0 


(m* = 0,l»2,...,M coef -l equations) 


Equation (12) comprises M co€ j- separate difference equations, each 

of which is written in terms of all the unknown coefficients A + and 

m 

the pressure pj at the nodes on S. These equations are combined 
with the finite element equations presented below to form a matrix 
system which yields all the unknowns at once. 

The continuity of velocity (Baumeister and Dahl, 1987, Eq. (21)) 
requires that at the interface 


1 dPa 1 dp 

£ » s . e dnls - 


where n is the outward normal. This relationship is key to the 
modal element method, as indicated below. 

BOUNDARY CONDITIONS 

The pressure field is scaled appropriately by introducing the 
normalized incident wave in Eq. (6). The far field radiation boundary 
condition is satisfied exactly by the assumed form of p® in Eq. (8). 
The pressure continuity across the domain interface is enforced by 
Eq. (12), and the velocity continuity across the interface by Eq. (13). 
No special handling of acoustic boundary conditions at media inter- 
faces within the finite element region is required because the variable 
property form of the wave equation has been employed. Any change 
in internal structure in the finite element region can be simulated 
merely by changing the properties of the elements in the finite ele- 
ment domain, as shown by Baumeister and Dahl (1989). 

FINITE ELEMENT SOLUTION 

The finite element domain, with total area A, is divided into M 
discrete triangular elements, A^ e \ e = 1,2,...,M, defined by N corner 
nodal points (x:,yj, i = ,2,...,N. The corner nodes for area A^ 
are denoted [ x i jYi ]» [x^,y^], and [x^*\y| e ^]. It is assumed that 
all material properties are constant in each element. 

The acoustic pressure is approximated by a linear combination 
of weight functions Wj(x,y): 



Figure 4.— Number of modes required to resolve the 
scattered acoustic pressure field around a cylinder of 
radius a subjected to a plane wave impingement, at 
various radii, for relative error e r = 10“®. 
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( 19 ) 


N 

p(*>y) = E w i(x.y)Pi = [w(x,y)]{p} , 

i=l 


(14) 


with [ ] representing a row vector and { } representing a column 
vector. The weights have the property that 


WifcpYj) = <*ij (Kronecker delta) , 


(15) 


so that the unknown nodal pressure vaJues are given by pj = p (x^yj). 

To determine {p}, apply the method of weighted residuals. In 
this method, the residual error of Eq. (1), 


R. = J J A Wj |v • Ivp + oi 2 /ipjdxdy 

(i * 1,2,3,...N one residual 

equation for each node i) 


(16) 


is set to 0 for each node i. Applying the divergence theorem to the 
first term of the integrand in Eq. (16) yields the weak formulation of 
Eq. (1): 




-Vp 

€ 


+ Wj^Vp 

dxdy - J 

/ 

WjlVp • n 


J b 

e 


ds = 0 . 


(i * 1,2,3,...N one global finite difference 

equation for each node i) (^) 


pW(x,y) = N| e, (x,y)p[ e) ♦ N< e) (x, y )p< e) + N< e) (x,y) P ' e) 
= £ Nj e) (x,y) P ; e) « [N<'>(x,y)j{p( e >} . 


Now implement the Galerkin method — let the local weight function 
Wj( e ) associated with node (xj^,y-^) equal N-^. The global 
shape function W. is identically zero for any element where node i 
does not appear. 

In the boundary integral terms in Eq. (18), it is reasonable to 
approximate the (continuous) normal derivative with its mean value 
over S^ViS. The key step is to apply the continuity of velocity 
(Eq. (13)), which introduces the eigenfunction coefficients, thus link- 
ing the analytic solution and the finite element solution on the inter- 
face. The term is transformed as follows: 


Is<*> ns ( W i e> ) V P W ' “jfc ■ Js ( .| fi! 

= 1 f N^ds 


ds 

1 e dn ; 


( 20 ) 


L f^cos(/? - 9) * i^sin(/J - 6 ) 

6 _l r de 


(<) Js('>nS N i ds 


where 6 and 0 are shown in the insert in Fig. 1. 

Substituting Eqs. (19) and (20) into Eq. (18) and employing 
the Galerkin approximation yield 


Equation (17) is a global, or node-oriented, formulation, in that it 
provides a difference equation for each node that can be used to 
determine {p}. 

From a practical standpoint, though, it is more convenient to 
consider a local, or element-oriented, formulation, which more easily 
provides the same difference equation for each node. To develop the 
local formulation, write each residual Rj as the sum of the element 
residuals: 


M M 

o = Ri=E R ! e> =E 

e— 1 e=l 


/ • - Vp (e) + wj e) wVpWjcU 


Js(') ns 


W. (e) IvpW ■ nlds 


(i = 1,2,3,...N one global finite difference 
equation for each node i) 


(18) 


where is the boundary of element A^. Note that for stand- 
ard weight functions, the line integral in Eq. (18) vanishes unless 
node i is on the boundary S. 

In order to evaluate the integrals in Eq. (18), it is necessary to 
represent p (x,y) locally. Let Nj^ e \ j = 1,2,3, be the local linear 
shape functions associated with each corner node (Segerlind, 1976, 
p. 29), so that 


M 

0 = Ri - E 

e=l 




J_VN? e) • V[N(*)| - wy e >N. (e, [N( e >] 

,(•) ‘ 


(p(')} dA 


^'_!*co3(/) 0) + i fl!?siil(/J 

dr r dO 


Js(*) ^ s{<)nsNi 

(i = 1,2,3,...N one global Finite difference 


( e >d S 


equation for each node 


( 21 ) 


For each value of i, Eq. (21) can be evaluated explicitly using the 
standard definition of Nj(x,y) for linear triangular elements 
(Segerlind, 1976, p. 29). 

Combining Eqs. (21) and (12) provides the final form of the 
matrix equation 


K {<!>} = {F}, 


where 


{^} T = |^i ’^2 , "‘ , ^M co<)f >Pl»P2»"**PN 


( 22 ) 


(23) 


F contains the incident plane wave terms present in Eq. (12) and in 
the derivatives of p a in Eq. (21). The matrix K has the following 
general form 
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A< 22 > 


i p(> 2 ) 

,_j 

: p( 22 ) . 


The submatrix A^ 12 ^ is a full M coef x M coef matrix composed of the 
coefficients of the terms in Eq. (12). The submatrix P^ 12 ^ is a 
sparse M coc j x N matrix composed of the coefficients of pj in 
Eq. (12). A^ 22 ) is an N x M coef matrix composed of the coefficients 

of A^ from the surface integral in Eq. (21). For each boundary node, 
there is a full row of terms in A^ 22 ), with a full row of zeros for 
interior nodes. P^ 22 ) is a sparse, highly banded N x N matrix 
composed of the coefficients of pj in Eq. (21). 

Equation (22) is solved by a frontal solver. The number of 
columns kept in core is roughly M coef + M fc0 /2 + bandwidth of 

p(22) p 

RESULTS AND COMPARISONS 


nificantly. For example, if ka = 100, then 110 modes are required 
for e r — 10 — \ but only 18 more modes are needed to improve the 
accuracy to e r = 10” 6 . Similar results hold off the surface, in both 
the near and far fields. Note that for the most computationally 
intensive case considered (ka = 100, e r — 10 6 ), the number of 
columns kept in core by the frontal solver is under 300. 

It might be expected that the number of required modes increases 
in the far field, but Fig. 4 indicates otherwise. As r increases, the 
amplitude of decreases, so the terms are smaller and the ser- 

ies converges faster. As seen in the figure, the decrease in the num- 
ber of terms is slight, and occurs even for small r/a values. The 
data in Fig. 4 is for e f = 10 3 , but is representative of all the other 
e r values. 

The following formulas give a rough indication of the number 
of modes needed for convergence: 


( 




r 

1 

M c = int 

7 + lo 6l0 I 

+ 

l.ie + o.o4iog 10 

1 

j 

(ka) 


In order to validate the method, several numerical experiments 
are presented. First, a description of the number of modes required 
for convergence in Eq. (8) is given, and a simple formula is provided. 
Example 1 shows the simple case of a plane wave passing through air. 
Examples 2 and 3 describe the problem of scattering from a hard rigid 
cylinder, using the modal element method and the modal ring method, 
respectively. Both numerical results agree closely with the exact 
solution for these problems. Finally, scattering from a hard circular 
cylinder with an absorbent coating is presented in Example 4. These 
examples give an indication of the value of the modal element method: 
the speed and accuracy of the calculations, the relatively low matrix 
storage costs, the wide frequency range available, and the flexibility 
to handle problems of complicated geometry and material property 
variation. 

For the circular cylinder examples to follow, the characteristic 
length is set equal to the radius of the circular cylinder. Conse- 

quently, the dimensionless radius a = 1 in all examples. 

It is of practical interest to determine the number of modes 
necessary for the series in Eq. (8) to converge. As more modes are 
needed, more finite element nodes on the interface are also needed, 
and hence more computing resources are required. To get a rough 
idea of the number of modes needed for a particular application, 
consider the expression for the scattered field generated by an inci- 
dent plane wave (as in Eq. (6)) striking a hard cylinder. This 
expression is given by the infinite series associated with Eq. (8), 
where (Bowman et ah, Eq. (2.38)) 


A T = 


-(2 - * m0 ) i m J m (ka) 


(24) 


ka < 15 , 




r ‘ 




1 

M^. = int 

10 + log I0 

l 

[e 

+ 

1.00 + 0.03log 10 

1 

(ka) 


ka > 15 . 


Equation (27) is valid for ka up to 100, and may be accurate for 
larger ka values. In many cases, the number of terms indicated by 
the formulas is slightly larger than the number of terms actually 
needed, but is nowhere smaller. The worst observed difference is for 
ka = 100, e r — 10 -6 , where 128 terms are required, but the formula 
gives 134. The equations are valid on the surface as well as in the 
far field. 

Figures 3 and 4 are not valid for a general problem. However, 
they represent a good starting point for estimating the number of 
modes required for convergence. 

In each of the following examples, a unit plane wave, incident 
from the left, strikes a cylinder oriented with its axis normal to the 
propagation direction. 

Example 1. - Propagation Through Air 

Consider the incident plane wave propagating through air and 
impinging on an imaginary cylinder (radius 1) composed of air. The 
modal element method verifies that there is no reflection from the 
cylinder and the pressure amplitude inside and outside the cylinder 
is 1. Figure 5 shows the pressure amplitude at the interface between 
the finite element and homogeneous domains. For a relatively coarse 
grid, the first few calculated modal coefficients (which should vanish) 
are: 


Clearly, the number of modes is determined by ka, the observation 
radius r, and the desired accuracy of the solution. For the data 
presented in Figs. 3 and 4, accuracy is measured relatively; that is, 
letting the partial sums of the series be denoted by S n , the series 
has converged if 


S n 

In Fig. 3, the number of modes needed for various e r values is 
compared with ka, on the surface of the cylinder (r = a). The data 
clearly shows that as e f decreases, the number of required modes 
increases, as expected. However, the increase is slow enough that 
only a few extra terms are needed to improve the accuracy sig- 


(25) 


A* = -0.33424 xlO -4 - 0.57905x 10~ 2 i 
A* = -0.13103x10' 2 - O.lSSOSxlO -4 i 
A* - 0.12062x10 5 + 0.45241xl0 _3 i 

Example 2. - Rigid Cylinder, Modal Element Metho d 

Consider the incident plane wave striking a solid rigid cylinder 

of dimensionless radius r = 1, with ka = 1. For all internal clc- 

19 - 

ments, the permittivity is set to e = 1 — 10 i to approximate a 
hard cylinder. The finite element domain extends to r = 1.5; in this 
case the polar region contains air, although it could be filled with an 
absorbing material, as in Example 4. 
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□ 


Exact analysis ABS (p) 
Numerical solution ABS (p) 


Clearly, the numerical and exact values of the modes are very close 
for those modes which make a significant contribution to the scat- 
tered field. 



Figure 5. — Polar plot of the total acoustic pressure 
around a cylinder composed of air (zero scattering, k 
= 0.628, r = 1 .0, 57 nodes, 96 elements, Ax = Ay = 
0.25 in rectangular region), interface at r = 1 . 


In Fig. 6, the pressure amplitude is plotted against azimuthal 
angle for ka = 1. Figure 6(a) shows the pressure on the cylinder’s 
surface (r = 1), while Fig. 6(b) shows the pressure at r= 1.5. Since 
these radii are both in the finite element domain, the finite element 
solutions (hollow boxes) are compared to the exact solution (solid line). 
Clearly, the method gives excellent results. 

In Fig. 7, the pressure amplitude is plotted again for ka = 1, 
but for larger radii in the homogeneous domain. In these plots, the 
numerical solution is generated from the calculated modal coefficients 
using Eqs. (5), (6), and (8). Plots for r = 1.5 (7a), r — 3.18 
(7b), r = 7.95 (7c), and r = 15.9 (7d) are shown. Again, there is 
excellent agreement between the numerical solutions (hollow boxes) 
and exact solutions (solid lines). 

The grid shown in Fig. 1 was used to generate the data pre- 
sented in Figs. 6 and 7. Recall that from the analysis summarized 
in Fig. 4, if a given set of modes accurately describes the pressure 
field on the surface of the solid cylinder, then the pressure is accu- 
rately described for all larger radii, since the required number of 
modes decreases in the far field. The first few calculated modal 
coefficients are: 


Numerical Solution 

A + = -0.23066E-Q0 - 0.42154E-0Oi 
A+ = — 0.67095E-00 - 0.256S3E-OOi 
A j = +0.14075E-01 - 0.166S8E-00i 
A 3 = +0.70539E-02 4 0.11636E-03i 
A+ = -0.49985E-06 4 0.1370lE-03i 
A+ = — 0.66943E-05 - 0.45628E-06i 
Ag = +0.13118E-07 - 0.21266E-06i 
A^ = — 0.55536E-07 - 0.29859E-07i 
Ag = +0.22121E-09 - 0.70338E-09i 


Exact Solution 

Aj = -0.24087E-00 - 0.4276lE-00i 
A j = -0.65615E-00 - 0.24537E-00i 
Aj = +0.1382SE-01 - 0.16569E-0Oi 
A 3 = +0.71092E-02 4 0.2527lE-04i 
A+ = — 0.11510E-07 + 0.15172E-03i 
Ag = — 0.19355E-05 - 0.1873lE-lli 
Ag = 40.15400E-15 - O.1037lE-O7i 
A+ = 40.98529E-10 4 0.48540E-20i 
Aj = — 0.98392E-25 4 0.44S60E-12i 


Exact analysis ABS (p) 

□ Numerical solution ABS (p) 




Figure 6. — Polar plot of the total acoustic pressure at 
node points around a solid cylinder subjected to a 
plane wave impingement (k = 1 .0, a = 1 .0, 1 865 nodes, 
3648 elements, Ax = Ay = 0.05 in rectangular region), 
interface at r = 1 .5. 
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Exact analysis ABS (p) 

□ Numerical solution ABS (p) 



(c) r = 7.95. (d) r = 15 . 9 . 

Figure 7. Polar plot of the total acoustic pressure around a solid cylinder subjected to a plane wave impingement 
reconstructed from modal-element method modal coefficients (k = 1.0, a = 1.0, 1865 nodes, 3648 elements, Ax = 
Ay = 0.05 in rectangular region), interface at r = 1 .5. 



□ 



0 


Example 3 - Rigid Cylinder, Modal Ring Method 

Consider a high frequency incident plane wave striking a solid 
rigid cylinder of dimensionless radius r = 1. The modal ring method, 
using a grid similar to that of Fig. 2, is used to generate the pressure 
amplitude plots in Fig. 8 for frequencies corresponding to ka — 50 
and ka = 100 at r = 1 and r = 7.9575. Choosing as an acceptable 
relative error e — 10 3 , the number of modes required to resolve 
the scattered field is M c — 67 for ka — 50 (Figs. 8(a) and (b)) and 
M c = 121 for ka = 100 (Figs. 8(c) and (d)). The total number of 
finite element nodes used is 24M C (12M c each on the interface and 
the outer ring). The excellent agreement between the numerical 
solutions (hollow squares) and the exact solutions (solid lines) 
clearly indicates that the modal ring method is suitable for high 
frequency scattering applications. 

As an additional consideration, when the modal ring method is 
applied to a low frequency scattering problem considered by Khan, 
Brown and Ahuja (their Fig. 3), the modal ring method gives more 
accurate results with a reduction of two orders of magnitude in the 
number of grid points. 

Example 4. - Coated Rigid Cylinder 

Consider the incident plane wave impinging on a cylinder of 
dimensionless radius r = 1 that has a rigid solid core of radius 0.9 
and an absorbing coating of thickness 0.1. The solid rigid core is 
approximated by setting e = 1 - 10 19 i while the coating has 
e = fi — 1 — 5i. The grid shown in Fig. 1 was used in this example. 
Figure 9 shows the calculated pressure amplitude at r = 3.18. The 
effect of the coating is apparent by comparing this figure with 
Fig. 7(b). The profiles have the same general shape but the coating 
induces slightly more backscatter (0 — 0), and slightly less scattering 



Figure 9. — Polar plot of the total acoustic pressure 
around a coated cylinder subjected to a plane wave 
impingement reconstructed from a modal-element 
method modal coefficients (k = 1 .0, a = 1 .0, r = 3.1 8, 
€ = 1 -i5, |x = 1 -i5, 1145 nodes, 2208 elements, Ax = 
Ay = 0.05 in rectangular region), interface at r = 1 . 


at 0 — 135, 225°. Because the modal element method is relatively 
inexpensive to run, it could be used to identify coatings that 
maximize the absorption of acoustic energy. 

CONCLUDING REMARKS 

The modal element method for acoustic scattering from a two- 
dimensional body is presented. The body may be acoustically hard 
(reflecting) or soft (absorbing). The acoustic pressure field is repre- 
sented by a finite element solution in the bounded finite element 
domain containing the body, and by an exact eigenfunction expan- 
sion in the unbounded exterior homogeneous domain. The two repre- 
sentations are coupled by the continuity of pressure and velocity across 
the interface between the two subdomains, and are calculated simul- 
taneously from a single matrix equation. The method is applicable to 
problems involving high or low frequency scattering. 

The analysis for two-dimensional scattering from solid and 
coated bodies is presented, and several simple numerical examples 
are discussed. The examples are provided for validation — the 
numerical results show excellent agreement with the corresponding 
exact solutions. More detailed numerical work is planned. 

For scattering from rigid bodies, the modal ring grid offers an 
attractive low storage approach to handling high frequency scatter- 
ing. The modal ring grid effectively reduces the two-dimensional 
problem to a one-dimensional problem similar to the boundary ele- 
ment method. The ring grid offers the added ad vantage of being 
able to handle rigid bodies with coatings. 
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